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Abstract 

We have performed extensive simulations on the Abelian Sandpile Model (ASM) on square lattice. 
We have estimated the probability of observation of many clusters. Some are in good agreement with 
previous analytical results, while some show discrepancies between simulation and analytical results. 

1 Introduction 

Abelian Sandpile Model (ASM) is a generalization of BTW model, proposed by Bak, Tang and Wisenfeld 
£Q for understanding the mechanism of Self-Organized Criticality (SOC). Although the physical relevance of 
ASM can be questioned, it is the simplest and most studied one among all other SOC systems. It is believed 
that, the model is an exactly solvable with all the features of the SOC systems 2 . There are very good 
reviews on this model, for example see 0E]- 

BTW sandpile model is defined as a simple cellular automata on the square lattice of size L x M. On 
each site i a variable hi is defined as its height, taking values {1,2,3,4}. System evolves in discrete time: 
at each step one site is picked randomly and its height is increased by unity. If its height becomes larger 
than h c — 4, this site is said to be unstable and it relaxes by toppling. The process of toppling is defined as 
follows: four grains leave the unstable site and each of its four neighbors gets one grain. So the number of 
grains is conserved during the toppling, except at the boundaries where some of the grains leave the system. 
If there is any unstable site remaining, it will also topple. 

The process of toppling could be stated in another way. If the site i becomes unstable, hj will be 
decreased by amount of Ay , that is hj — > hj — Ay where 

4, i=j; 

Aij = ^ —1, i,j are neighbors; (1) 
0, otherwise. 

After a while the system reaches a steady state, in which it shows SOC. It has been shown that in this 
state, the number of different configurations the system accepts is det A and the probability of all of them 
are the same. These configurations are named recurrent configurations in contrast with the transient ones 
which can only appear in the first steps of evolution, where the system has not yet reached the steady state. 

The question is, what is the probability of finding one given cluster in SOC mode? It is possible to show 
that some sub-configurations are disallowed in SOC mode (e.g. clusters in fig-H-a). There is an easy way 
to test whether a given cluster is forbidden or allowed, so the probability of some clusters are exactly equal 
to zero. Dhar and Majumdar calculated the probability of finding some allowed clusters using a beautiful 
idea [5]- This method is appropriate only for weakly allowed clusters. These clusters have the property that 
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decreasing the height of any of their sites, makes them forbidden (such as So in fig. ^b and all clusters in 
fig. ^-c and fig.njd). The first cluster calculated analytically in [S] was Sq: 

P(S ) = A - A = 0.0736362... (2) 

A table of some other analytical results can be found in 0. It is obvious that 5io and Sn in ref. |U] are not 
weakly allowed clusters and one can not use this method for them. 

Yet there are many other clusters that are not included in this set of clusters, so their probabilities are 
still undetermined. It seems that some generalization is required for exact calculation of probability of these 
clusters. Some attempts have been made in this direction, for example in jB] all height probabilities have 
been calculated. 




Figure 1: Some different types of clusters: (a) Forbidden Sub-Configurations, (b) Allowed Clusters (So is also 
Weakly Allowed Cluster), (c) S$, Sg have exactly the same analytical values but different simulation results, (d) 
Some of weakly allowed clusters which are in disagreement with analytical results. 

Several works have been done based on analytical approach in [5] such as calculation of correlation 
functions 0, probability of clusters and correlation functions in boundaries |S] and scaling fields 
Nevertheless, we focus only on the first problem: "the probability of weakly allowed clusters" because it 
contains all foundations of analytical approach throughout the literature. 

2 Simulation results 

The first simulations were done on small lattices and for simple clusters such as Sq. Extensive simulations 
for the lattice of size 672 x 672 were undertaken by Manna in [?], who found -P(Su) = 0.0736 ± 0.003. 
Grassberger and Manna also performed simulations on some larger lattices, but not with sufficiently high 
statistics 

Although these results are in agreement with analytical result in eq. (0, it seems more clusters should 
be tested with various sizes and forms. It is also needed to perform more precise simulation using larger 
lattices and more samples. 

We first started with a lattice of size 500 x 500 and extended it to a 2000 x 2000 one and measured the 
probability of clusters introduced in the table of paper 9 for 2,000,000 samples, sampling is started when 
the system reaches SOC using the burning test |12| . Each new sample is obtained by adding a grain of sand 
and allowing the system to relax. This process is repeated for 2,000,000 times and we find a distribution 
function for probability of appearance of any given cluster. The result is typically a Guassian distribution 
(see fig. 01. 

Statistical analysis shows that for one of the clusters (5*8) there is a deviation from analytical result (fig. 0- 
a) 1 . To confirm that the deviation does exist and is not a simulation error, we have done the same simulation 

x In lattice of size 500 X 500 , the deviation was not so clear, though one could guess that there is something wrong. 
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for the cluster Sg , which topographically is the same as Sg . Also probability of Sg is analytically the same 
as that of 5*8 and simulation confirms its analytical result without any remarkable deviation as fig- EJb. 

ft is also noticeable that the quantity studied by analytical approach, is an average over all recurrent 
configurations. Standard deviations in fig. J2J-a and fig. 0-b are the error in each individual measurement of 
the sample, so the more appropriate quantity to stand for error bar is the standard error of the mean (<r m ). 
fn our simulations it was generally much less than the error of precision concerning each individual sample 
which is equal to jj where N is the number of all attempts to find the specific cluster in each sample. This 
quantity, therefore, is the best one for estimating the error. 

The simulation results for the clusters Sg and Sg were obtained : 



P 2 ooox2ooo(5 , 8 ) = 0.000162467 ± 0.000000031 a = 0.000003072 a m = 0.0000000022 
^2000x2000(59) = 0.000173098 ± 0.000000031 a = 0.000002558 a m = 0.0000000019 

to be compared with the analytical result : 

PAnalytic(5 8 ) = P A nalytic (Sg ) = 0.000173106... 



(3) 
(4) 



(5) 



We repeated the simulation on large lattices of size 500 x 500, 1000 x 1000 and 2000 x 2000 for all clusters 
of size < 6. We found that there are more clusters in disagreement with the analytical results (e.g. clusters 
in fig. ffl-d), though they are in a minority. 

Repetition of simulation confirms our results and shows that there is a noticeable disagreement between 
analytical result and simulation. The origin of this disagreement remains to be understood. 





a) Probabiliy of S 8 



(b) Probabiliy of S 9 



Figure 2: Distribution of probability of clusters over 2,000,000 samples on lattice of size 2000 x 2000: (a) P(Ss) is 
in disagreement with analytical value (thick bar), (b) P(Sg) is in good agreement with analytical value. 



3 The source of deviation 

Either simulation circumstances are not appropriate or something in analytical approach is not sufficiently 
reliable, fn the former, some possible questions are imaginable and we discuss them to rule out this possibility. 

One proposal is to do the simulation with more samples, because the number of recurrent configurations 
is much more than something like 2,000,000 (det A ~ (3.21) ixM ) and so, our samples cover only a small 
part of configuration space of system. To rule out such a possibility we repeated simulation for a variety 
of samplings from 100,000 to 10,000,000 samples and no meaningful change was observed in the amount of 
deviation, fn fact, the ergodicity of model makes our work reliable. Furthermore for each sampling, the 
initial recurrent configuration is generated randomly but the behavior observed in all of them is the same. 
Therefore this deviation cannot specific to some particular parts of configuration space. 

Other critique of our claim may be the difference between the definition of probability of cluster in 
analytical approach and in our simulations. In analytical method, probability is defined in terms of the 
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number of configurations in which the specific cluster appears in a given place, but in our case it is defined in 
terms of the number of appearance of the specific cluster in a given configuration. Although this difference 
really exists, the average of probability over all places in the former should be equal to the average of 
probability over all recurrent configurations in the latter. So this can not address our question about the 
source of deviation. While this argument seems to be enough, we measured probability of clusters Sg and 
Sg directly with the first definition in different areas. In all cases we found the same difference as before 
between P(Sg) and P{S 9 ). 

The other possible source of error could be finite size of the lattice in simulations, as the analytical 
values have been calculated exactly for infinite lattice. It is clear that for very small lattices, the calculated 
quantities do differ noticeably with results of infinite lattice. In larger latices, finite size effect may also cause 
a little deviation from infinite case, but this can not be the source of observed deviation. It seems that the 
finite size effect is much smaller than the observed deviation in our simulation, as it has not appeared in 
majority of clusters in our simulation circumstances. 

Another point about this critique is the fact that the similarity of Sg and Sg does not concern only 
infinite lattice. According to the analytical approach, the probability of appearance of clusters 5*8 and Sg 
are equal in every lattice and this claim of analytical approach is in considerable disagreement with our 
simulation results. 

Using simulation for lattices with different sizes from 5 x 5 to 2000 x 2000 shows that the difference 
between P(Sg) and -P(Sg) not only does not decrease with size, but also increases; whereas for lattice of size 
5x5, the difference between the two values is not statistically significant. But the problem arises quickly 
and becomes more serious for larger lattices specially because the error bar narrows down. 

Finally, we may now focus on the analytical results. There does not seem to be a calculational problem. 
Therefore we suspect that a fundamental problem is involved which is however obscure to us. 
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